Optimal Bounds on the Trapping Constant and Permeability of 

Porous Media 

S. TorquatcQ and D. C. Pham 

Department of Chemistry and Princeton Materials Institute, 
Princeton University, Princeton, NJ 08544 
(Dated: February 2, 2008) 

Abstract 

We derive exact expressions for so-called "void" bounds on the trapping constant 7 and fluid 
permeability k for coated-spheres and coated-cylinders models of porous media. We find that in 
some cases the bounds are optimal, i.e., the void bounds coincide with the corresponding exact 
solutions of 7 and k for these coated-inclusions models. In these instances, exact expressions are 
obtained for the relevant length scale that arises in the void bounds, which depends on a two-point 
correlation function that characterizes the porous medium. In contrast to bounds on the effective 
conductivity and elastic moduli of composite media, this is the first time that model microstructures 
have been found that exactly realize bounds on either the trapping constant or fluid permeability. 
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The study of the effective properties of heterogeneous materials, such as composite and 



porous mediajhas a rich history 
questions Q, E, B, S fl 2, 3 



03 



and is a continuing source of theoretical' 



y challenging 



exact 



Except for a few special microstructures |6, 
predictions of the effective properties are not possible because they depend on an infinite 
set of statistical correlations that characterize the microstructure [8j]. Thus, apart from such 
exact solutions, rigorous estimates of effective properties must take the form of inequalities, 



i.e., upper and lower bounds 



QBQ r 



lllj . Bounds are useful because often one of 



the bounds can provide a useful estimate of the property even when the reciprocal bound 

n 

diverges from it Moreover, it is highly desirable to determine optimal bounds and the 
microstructures that attain them. The best known bounds in the cases of the effective 
conductivity and bulk modulus of two-phase media are the Hashin-Shtrikman bounds [lGj. 
These are optimal bounds, given the phase volume fractions, because they are realizable by, 
among other geometries, certain coated-spheres and coated-cylinders assemblages in three 
and two dimensions, respectively. 

Two important effective properties of fluid-saturated porous media that have been ex- 
pensively studied are the trappingconstant 



Bounds on 7 




and scalar fluid permeability k 



0, 22,0, 24,0 



have been 



derived and computed. However, to date, microstructures that exactly realize (or attain) 
any of these bounds have yet to be identified. Torquato 8] has observed that the difficulty 
in identifying optimal microstructures for these classes of problems lies in the fact that 7 
and k (unlike the conductivity and elastic moduli) are length-scale dependent properties and 



known bounds on them depend nontrivially on the specific forms of two-point anc 
order correlation functions. For example, the so-called void bounds on 7 and k j^, 



ligher 



21 



depend on the two-point correlation function 5*2 (r) (defined below) and have been evalu- 



ated for various particle models for the trapping constant |8 



21 



and fluid permeability 



In this Letter, we exactly evaluate the void bounds on the trapping constant 7 and 
fluid permeability k for the coated-spheres and coated-cylinders models of porous media. 
Interestingly, we show that in some cases the void bounds are optimal because they coincide 
with the corresponding exact solutions of 7 and k for these particular coated-inclusions 
porous-media models. In these cases, we obtain exact expressions for the relevant length 
scale that arises in the void bounds, which depends on S2 (r). 
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Each realization of the porous medium occupies a region of d- dimensional space V of 
volume V that is partitioned into two disjoint regions: a pore space (phase) Vp of porosity 
(ftp and a solid space (phase) Vs of volume fraction <p s = 1 — <$>p. Let dV denote the surface 
of interface between Vp and Vs- The pore-space indicator function X^(x) is given by 

Z«(x) = ( 1 ' XmVp (1) 
I 0, otherwise 

The indicator function .M(x) for the interface is defined as 

M{x)=\ VJ (p) (x) | • (2) 

For statistically homogeneous media, the ensemble averages of the indicator functions ([TJ 
and (J2J) are respectively equal to the phase volume fraction (ftp and the specific surface s 
(interfacial area per unit volume), i.e., 

0P=(j( p >(x)>, s = {M{x)) } (3) 

where angular brackets denote an ensemble average. 

Before evaluating the bounds, we first define the effective properties 7 and k, and the 
coated-spheres model. First, consider the steady-state trapping problem ^]. The reactant 
diffuses in the pore space Vp (i.e. trap-free region) with scalar diffusion coefficient D but is 
instantly absorbed when it makes contact with the interface between Vp and the trap region 
Vs- At steady-state, the rate of production of the reactant G (per unit trap- free volume) 
is exactly compensated by its removal by the traps. Two-scale homogenization theory 20| 
enables one to show that the trapping constant 7 for a statistically homogeneous and ergodic 
medium obeys the first-order rate equation 

G = jDC . (4) 

Here C is the average concentration field and 

T l = («) = J im 77 f m(x)c?x , (5) 

V 

where w(x) is the scaled concentration field that solves the boundary- value problem 

V 2 w(x) = -1 , x G Vp , (6) 
u(x) = , x G dV . (7) 



It follows that the trapping constant 7 for any d has dimensions of the inverse of length 
squared |8|. 

n 

Rubinstein and Torquato [20J have formulated a variational principle in terms of the trial 
function u(x) that enables one to obtain the following lower bound on 7 for ergodic media: 

7 - (Vv(x) • Vw(x)X^)(x))' (8) 
where t>(x) is required to satisfy the Poisson equation 

V\;(x) = -1 , x G V P . (9) 



Elsewhere Torquato and Rubinstein |21| constructed what they referred to as the void lower 
bound in three dimensions by using a specific trial field. The generalization of this trial field 



to any space dimension d > 2 



is given by 



«(x) = -J- / </(x- y)[x (p) (y) - <M^y, (10) 

where 

9{r) = \ 1 1 d>3 ( ] 

is the d- dimensional Green's function for the Laplace operator [8], Q(d) = 2-K d l 2 /T(d/2) is 
the total solid angle contained in a <i-dimensional sphere, 05 is volume fraction of the trap 
phase, and r = |r|. Substitution of trial field (jTU|l into the variational principle (jHJ) yields 
the two-point void lower bound on 7 for general statistically homogeneous and isotropic 
d-dimensional porous media |8] as 

7>f, (12) 

Ip 

where ip is a pore length scale defined by 

-/(TWO - <p 2 P ]r In rdr, d = 2 



f2 



(13) 



d > 3, 



^2jJo°[^(r) - ^,]rdr, 
and ^2 (r) is the two-point correlation function defined by 

S 2 {r) = (T^(x)Z^(x + r)) . (14) 
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The function S2 (r) can also be interpreted as being the probability of finding two points 
separated by the displacement vector r in the pore space 

Using homogenization theory, Rubinstein and Torquato |25( derived the conditions under 
which the slow flow of an incompressible viscous fluid through macroscopically anisotropic 
random porous medium is described by Darcy's law 

U = --k-V Po , (15) 

where U is the average fluid velocity, Vpo is the applied pressure gradient, fi is the dynamic 
viscosity, and k is the symmetric fluid permeability tensor. In particular, for the special 
case of macroscropically isotropic media, the scalar fluid permeability k = Tr(k)/d (where 
Tr denotes the trace operation) is given by 

k = (w ■ e) = lim — \ w(x)<ix , (16) 

where e is a unit vector, and w and it are, respectively, a scaled velocity and scaled pressure 
that satisfy the scaled Stokes equations 

V 2 w = Vtt - e in V P , (17) 
V ■ w = in V P , (18) 
w = on dV . (19) 

It follows that the permeability k for any d has dimensions of length squared 



Prager 



22j was the first to derive a two-point "void" upper bound on the permeability 



using a variational principle. Subsequently, Berryman and Milton 23( corrected a nor 



malization constraint in the Prager variational principle using a volume-average approach. 
Rubinstein and Torquato 3| formulated new upper and lower bound variational principles 
employing an ensemble-average approach and also derived the void upper bound. 

For our purposes, the Rubinstein- Torquato variational principle for the upper bound is 
the most natural. It states that for ergodic media the trial function q(x) enables one to 
obtain the following upper bound on A;: 

K(q(x):Vq(x)l( p )(x)), (20) 

where q(x) is required to satisfy the momentum equation 

VxV 2 (q + e)=0, xGVp (21) 



Rubinstein and Torquato 25] constructed the void upper bound in three dimensions by 
using a specific trial field. The generalization of this trial field to any dimension d > 3 is 



given by 



q(x) = J_y *(x- y) • e[J^(y) - <p P ]dy, (22) 



where 

* (r)= (^-3MV-- |I + H - d - 3 ' (23) 

is the (i-dimensional Green's function St (second order tensor) associated with the velocity 
for Stokes flow, n = r/r, and cj>s is the volume fraction of the obstacles. Substitution of trial 
field (j22j) into the variational principle (J2Uj) yields the two-point void upper bound on k for 
general statistically homogeneous and isotropic d- dimensional porous media [§! as 

where lp is the length scale defined by ([13)1 . which is precisely the same as the one that 
arises in the void lower bound on the trapping constant 7 for d > 3 [29^ . 

The coated-spheres model jlfl] consists of composite spheres that are composed of a 
spherical core of phase 2 (inclusion) and radius Ri, surrounded by a concentric shell of phase 
1 (matrix) and outer radius Rm- The ratio (Ri/RmY is fixed and equal to the inclusion 
volume fraction 2 in space dimension d. The composite spheres fill all space, implying that 
there is a distribution in their sizes ranging to the infinitesimally small (see Fig. Q). The 
inclusion phase is always disconnected and the matrix phase is always connected (except at 
the trivial point 02 = 1)- 

The coated-spheres model places restrictions on the size distribution of the composite 
spheres. Without loss of generality, we will assume that the composite spheres possess an 
infinite number of discrete sizes. Let pk be the number density (number of particles per 
unit volume) of the kth type of composite sphere of radius Ru k and let Ri k denote the 
corresponding radius of the inclusion. Moreover, we know that fraction of space covered by 
the composite spheres, denoted by is unity, and therefore we have the following condition 
on the size distribution: 

00 

$ = £/W#mJ = 1, (25) 

k=l 
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FIG. 1: Schematic of the coated-spheres model microstructure. 



where Vi(r) is the <i-dimensional volume of a single sphere of radius r given by 

Tt d/2 

Vl ( r ) = — —r d . (26) 

v ; r(l + d/2) v ; 

For the volume fraction $ to remain bounded (i.e., for the sum (|23|) to converge), PkR d M k 
must also remain bounded for all k. Thus, the number density pk must diverge to infinity 
as RM k approaches zero and the specific surface s must also diverge, since PkR d u k diverges 
as Ru k approaches zero. Note that volume fraction 02 of the inclusion phase is given by 

oo 7r^/ 2 00 

It is convenient to introduce the following nth moment of Rf. 

1 oo 

<*?> = -£fltf?*> (28) 
P k=i 

where p is a characteristic density (e.g., inverse of the volume of the largest composite sphere) 
and n is any integer n > 3. The inclusion volume fraction 02 can now be reexpressed as 

^ = "F(itW^ (29) 

We first evaluate the void lower bound on 7 for the three-dimensional coated-spheres 
model. To begin, we take the connected matrix phase Vi to be the traps and the disconnected 
inclusion phase V2 to be the pore space. Therefore, the porosity is given by <pp = 02- Using 



7 



the void trial field (fTU|) for v(x), we can obtain from (JHJ) the following lower bound on 7: 



7> 



1 f 

lim — / Vv ■ VvJ (2) (xWx 
^00 V J 



v 

v 



-1 



(30) 



where we have equated ensemble averages with volume averages via the ergodic hypothesis. 
The key v olum e „ ca n be ev.ea.ed Mowing P ham 0, and SU cn detaUs b e 
given elsewhere j^J. We find that the void lower bound is exactly given by 

„ 15(-R?) ,,,, 

"-MW)- (31) 

Interestingly, by comparing this result to the general expression for the void upper bound 
(112}) . which is given in terms of the two-point correlation function S2 (r), we see that the 
square of the pore length scale £p for d = 3 is exactly given by 

for the coated-spheres model. 

Now we show that bound ()31jl coincides with the exact solution for this particular coated- 
spheres model. Specifically, the exact solution of the boundary-value problem (JEJ) and (J7J) 
for diffusion inside a spherical inclusion Sj of radius i?/ is given by 8] u = (Rj — r 2 )/6 for 
< r < Rj. Hence, using definition (jSj), we find that 7, for nonoverlapping sphere models 
with a general size distribution (not just the coated-spheres model) is exactly given by 



7=<«> -1 



1 f R i 1 

,! im 77 J2 / «(^/ -r 2 Hiir 2 dr 



-1 



(33) 

0p(i?f)' (33) 



Thus, we see that the void lower bound (|3T|) coincides with the exact solution (|3"3"|) for 
the coated-spheres model, and hence the bound is exactly realizable when the inclusions are 
taken to be the pore phase. This may immediately lead one to conclude that the void bound 
is optimal among all microstructures, but such a statement cannot be made unless one 
attaches special conditions. Recall that unlike the effective conductivity or effective elastic 
moduli, the trapping constant as well as the fluid permeability are length-scale dependent 
quantities. Thus, any statement about optimality must fix not only the porosity but the 
relevant length scales. The correct statement is the following: The void bound is optimal 
among all microstructures that share the same porosity (ftp and pore length scale ip as the 
coated-spheres model [cf. Eq. (JH2J)]- One can always adjust the pore length scale (j3*2*|) of 
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the coated-spheres model at some porosity (ftp to be equal to i P for any microstructure with 
same porosity. 

As noted above, relation (|33|) applies to diffusion within nonoverlapping spheres with 
a general size distribution. Accordingly, let us define another squared length scale Lip = 

) / (Rj) for such a general nonoverlapping sphere model. In what follows, superscripts 
g and c are appended to quantities associated with the general sphere model and coated- 
spheres model, respectively. The use of expressions ffl"2j) . (J3"T|) . and (J3~3"j) reveals the following 
interrelations between these two models: at fixed 0p, if Lf = Lf, then 7^) = 7^ and 
if > if, and if if = if, then 7^ > and L ( f < hf . 

Next we take the connected matrix phase V\ to be the pore space and the disconnected 
inclusion phase V2 to be the traps so that (ftp = 4>i- Employing the trial field (jTUjl for v (x), 
we obtain the void lower bound for the model as 

->{l + UT-<Ps-Ul)-\ (34) 
7s o 

where 7 S = 3<fis(Rp)/ (Rf) and (f>s is the volume fraction of the traps. The exact solution in 
this case is not known. 

The procedure above can be repeated the two-dimensional coated-cylinders model. For 
diffusion inside circular inclusions, we obtain the void lower bound as 

which coincides with the exact result, and thus is an optimal bound in the sense described 
above. Comparison of this result to the general relation for the void upper bound (|12|) 
yields the following exact expression for the square of the pore length scale ip for the 
coated-cylinders model: 



e p = - jf> 2 (r) - <Pl]rlnrdr = *f^. (36) 



lo ™ 8 (R*y 

For diffusion exterior to the circular inclusions, we obtain the void lower bound as 

^>(-ln0 S+2^-^!)- 1 (37) 
7s 2 2 

where ls = 4(Rf)<p s /(R*). 

Consider fluid flow along (inside or outside) bundles of parallel cylindrical circular tubes 
corresponding to the coated cylinders model. The velocity field reduces to an axial compo- 
nent only, and the Stokes equation reduces to a simple Poisson equation identical to that of 



the 2D trapping problem. Hence, we have exactly the same solution for the axial component 
of velocity as for the concentration field in the trapping problem, leading to the exact result 
that k = 7 _1 [8]. Exploiting this observation and using the previous results, we simply 
summarize the appropriate results below for k. 

In particular, for axial flow inside the cylindrical tubes (Poiseuille flow), we obtain the 
void upper bound as 

k - ^RfT ' (38) 



which coincides with the exact result |8j, and thus is an optimal bound. Comparison of this 
result to the general relation for the upper bound (}2"4*j) also yields the exact formula (JHUj) for 
the square of the pore length scale ip for the coated-cylinders model in the transverse plane. 
A well-known empirical estimate for k is the Kozeny-Carmen relation (^J 

* = 4> ( 39 ) 
cs z 

where c is an adjustable parameter and s is the specific surface defined by 0. However 
for the coated-spheres or coated-cylinders model with the inclusions of all sizes down to 
infinitesimally small, we saw earlier that s diverges to infinity and therefore the Kozeny- 
Carmen relation incorrectly predicts a vanishing permeability. This serves to illustrate the 
well-established fact that the permeability cannot generally be represented by a simple length 
scale, such as the specific surface P- fl 

For flow exterior to the cylindrical tubes, we obtain the void upper bound as 

_ < _i n 5 _ 3 + 2 s - (4Q) 

?5\ KAM /f?3 



where k s = (Rj) /(Acp s (Rj)). For flow exterior to spherical obstacles, we exploit the fact that 
the void upper bound ()24|) on k is trivially related to the void lower bound (|12|) . Thus, we 
deduce the upper bound on k for flow exterior to spherical inclusions in the coated-spheres 
model from the corresponding bound (fH4"j) on the trapping constant: 

k/k s < 1 + ^cj>f -cPs- ^cf s , (41) 

where k s = 2{R]) / '(90 s (i2?». 

To summarize, the two-point void lower bound ()12j) on the trapping constant 7 and the 
void upper bound (|24*j) on the fluid permeability k both generally depend on the pore length 
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scale £p, defined by ()13jh which involves an integral over the two-point correlation function 
5*2 (r) that characterizes the porous medium. We have derived exact expressions for the 
void bounds on 7 and k for certain coated-spheres and coated-cylinders models of porous 
media. For diffusion inside the spherical (d = 3) and cylindrical inclusions (d = 2), the void 
lower bound on 7 was shown to be exact. Similarly, for axial flow inside the cylinders of 
the coated-cylinders model, the void upper bound on k was demonstrated to be exact. In 
these instances, the void bounds are optimal among all microstructures that share the space 
porosity 0p and pore length scale lp as the coated-spheres model. In contrast to bounds on 
the effective conductivity and elastic moduli of composite media, this is the first time that 
model microstructures have been found that exactly realize bounds on either the trapping 
constant or fluid permeability. For cases of diffusion and flow exterior to the spheres and 
cylinders in the coated-inclusions model of porous media, exact results are not available, 
but we still obtained simple analytical expressions for the void bounds on 7 and k. In 
future studies, it will be of interest to investigate what are the optimal microstructures that 
correspond to the improved two-point "interfacial-surface" bounds on both 7 and k, which 
also incorporate surface correlation functions 
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